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Abstract 

Shell Model Monte Carlo (SMMC) techniques are used to calculate two- 
neutrino double beta decay matrix elements. We validate the approach 
against direct diagonalization for ^'^Ca in the complete p/-shell using the 
KB3 interaction. The method is then applied to the decay of ''^Ge in the 
(0/5/21 059/2) model space using a newly calculated realistic interaction. 
Our result for the matrix element is 0.13 ± 0.05 MeV~^, in agreement with 
the experimental value. 

PACS numbers: 21.60.Cs, 21.60.Ka, 27.40.+Z, 21.10.Ma 
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The double beta (PP) decay of a nucleus is a rare second order weak process [jl],0]. The as 
yet unobserved neutrinoless mode is of fundamental interest, as it would signal a neutrino 
mass, lepton number non-conservation, or admixtures of right handed weak currents. In 
contrast, the existence of the 2z/ mode has been firmly established (see the review in Ref 
0). The ability to accurately describe this latter process is an important element in the 
interpretation of limits on neutrinoless decays. Unfortunately, it seems that 2z/ matrix 
elements are highly suppressed and so depend sensitively on small, poorly determined parts 
of the nuclear wavefunctions. 

Most recent calculations of 2i^/3/3 matrix elements for nuclei heavier than ^^Ca rely on 
the quasi-particle random phase approximation (QRPA) 0. While this approach is com- 
putationally simple and includes many features of the two-body interaction known to be 
relevant for (3(3 decay, the calculated matrix elements are uncertain because of their great 
sensitivity to the J = 1~^,T = particle-particle interaction The interacting shell model 
offers a more microscopic approach to the problem. Complete Ohu) shell model calculations 
^ not only recover more quenching of Gamow- Teller (GT) transitions than QRPA calcula- 
tions, but also are in agreement with observations (after the universal renormalization of qa 
to 1.0). However, computational limitations have restricted shell model calculations of the 
2u(3(3 decay matrix element to ^^Ca, the lightest of all (3(3 candidates. 

In this Letter we show how SMMC methods can be used to calculate 2v(3(3 decay matrix 
elements. We first calculate the decay of '^^Ca in the complete p/-shell and validate our 
method against direct diagonalization. We then present results for the decay of ^^Ge, one of 
the few nuclei where the 2v(3(3 decay has been precisely measured and where the best limits 
on the decay have been obtained 0. Our calculation in the (0/5/2,1^,0(79/2) orbitals, 
which is impractical using traditional shell model methods, is the first for 2v(3(3 decay in 
such a large model space. 

The 21/(3(3 matrix element between the 0^ ground states of the initial and final even-even 
nuclei is given by 0, 



2 



,,2._v (o;iG|i+)-(i;t.|G|o+) 

Here, |0j^) (|0j)) is the ground state of the initial (final) nucleus with energy (-^/)! I^m) 
are the 1^ states of the intermediate odd-odd nucleus with energies Em', and G is the GT 
operator, J^i^^i^r^ where ct/ is the Pauli spin operator for nucleon / and rf is the isospin 
lowering operator that changes a neutron into a proton. 

Previous shell model calculations for nuclei heavier than ^^Ca ( |l|] and references therein) 
have invoked the so-called closure approximation, where the matrix element is written as 
M'^" = Mc/E, with Mc = (Oj |G ■ G|Oj^) and E an average energy denominator. As there 
is no prescription for choosing the average energy denominator (and even the closure ma- 
trix element is usually calculated in a severely truncated basis), the uncertainty in this 
approximation is difficult to estimate. 

To calculate the exact 2i>(3f3 matrix element, Eq. (1), we consider the function 

„ Tr[e-(^— -')^Gt ■ Gte-^^Ge--'^ ■ G] 
^(-^ = (2) 

where H is the many-body Hamiltonian and the trace is over all states of the initial nucleus. 
The quantities {(3 — t — r') and r play the role of the inverse temperature in the parent and 
daughter nucleus respectively. A spectral expansion of shows that large values of these 
parameters guarantee cooling to the parent and daughter ground states. In these limits, we 
note that 0(t, r' = 0) approaches e~'^'^\Mc\'^ , where Q = E^ ~ Ej is the energy release, so 
that a calculation of 0(r, 0) leads directly to the closure matrix element. If we then define 

r/(T,r)= ^rfr'0(r,r')e-^'«/^ (3) 

and 

M-(T r) = (4) 

it is easy to see that in the limit of large r, (/3— r— r'), and T, M^'^(T, r) becomes independent 
of these parameters and is equal to the matrix element in Eq. (1). 



We use SMMC methods to calculate 0(r, r'), and hence M^'^. These techniques scale 
more gently than direct diagonalization with the number of valence nucleons and single 
particle orbits and so allow calculations larger then possible otherwise. They are based 
on the discretization of the many-body propagator, e~^^ , into a finite number of "time" 
slices, Nf, each of duration A/5 = P/Nf. At each time slice the-many body propagator is 
linearized via the Hubbard-Stratonovich transformation [§] ; observables are then calculated 
as expectation values in the canonical ensemble of nuclear states. 

To circumvent the "sign problem" encountered in the SMMC calculations with realis- 
tic interactions, we use the extrapolation procedure outlined in 0. One defines a set of 
Hamiltonians H{g,x) = (1 ~ ^^—^)Hc + qHb such that H{g = 1,%) = is the physical 
Hamiltonian and Hq^b are the "good" and "bad" parts of the Hamiltonian, respectively. 
For (? < 0, H{g^ x) is free of the sign problem; the matrix elements are therefore calculated 
for several values < and extrapolated io g = 1. The value of x is chosen to make the 
linear gf-extrapolation as smooth as possible. 

To validate our method, we calculated the matrix elements for ^^Ca in the complete pf- 
shell with the KB3 interaction for six equally spaced g values between -1.0 and 0.0 using 



X = 4 and extrapolated to the physical result ai g = 1. Each calculation involved 2500- 
3500 Monte Carlo samples and was performed at /3 = 2 MeV^^ with Nt = 48. The direct 
diagonalization calculations for ^^Ca with which we compare our results were performed 
using an implementation of the Lanczos algorithm [jll|. We calculated both the closure and 
exact matrix elements for the same Hamiltonians, H{g,x), a^s used in the SMMC. 

We found the slope of ln[0(r, 0)] to be in good agreement with that expected from the 
difference of the energies for ^^Ti and ^^Ca (Fig. 1) and extracted \Mc\ from the intercept. 
The SMMC closure matrix elements for g < are in very good agreement with the di- 
rect diagonalization results (Fig. 1) indicating that our temperatures are sufficiently low 
to correctly calculate the closure matrix element from the ground state of ^^Ca to ^^Ti. 
However, the direct diagonalization calculations show a small curvature near g = 1.0 that 
the extrapolation cannot reproduce. Our linear extrapolation of the closure matrix element. 



which takes place over almost a factor of 20, therefore underestimates the physical {g = 1.0) 
calculation. We obtain —0.21 ± 0.29 for the closure matrix element to be compared with 
the direct diagonalization result of 0.29. As the natural scale for Mc is given by the sum 
rule 10] as ~ 21, we may conclude that the SMMC successfully reproduces the shell model 
suppression of a factor of 70. 

The calculation of the function 0(r, r') was performed for r = 0.5 MeV^^ and for thirteen 
r' values spaced equally between 0.0 and 0.5 MeV~^. This combination of parameters was 
checked to give converged results for the matrix element. We then calculated ri{T, r) for 
thirteen values of T < 0.5 MeV~^; the upper limit of T is sufficiently large for the integral 
in Eq. (3) to converge. From these, we obtained M'^'^{T,t) (Eq. (6)) as shown in Fig. 2 



for some representative values of g [jT2[. In Fig. 2, we show the good agreement between 
the SMMC matrix elements and direct diagonalization for g < 0. Even though the value 
of X = 4.0 was chosen to make the linear extrapolation as smooth as possible, the direct 
diagonalization results still have a small curvature. For the exact 2z/ matrix element we 
obtain an extrapolated value of 0.15 ± 0.07 MeV~^ whereas the calculation of Caurier et al. 
]r3| (including the erratum in []TU[) gives 0.08 MeV^^. There is thus agreement within the 
uncertainty. 

We now apply the SMMC method to a heavier nucleus, where direct diagonalization is 
not possible. In particular, we calculate the 2v matrix element for ^^Ge using an effective 
interaction based on the Paris potential in the (0/5/2,1^,0(79/2) orbitals, with the single 
particle energies taken from the levels of ^''Ni relative to the ^^Ni core This interaction 



has been constructed using a G-matrix folded- diagram method, in close analogy with the 
calculations carried out by Shurpin et al. [|15| and by Dean et al. [Q. The model space 



comprises some 10^ configurations, so that our SMMC calculation is significantly larger 
than previous shell model treatments of ^^Ge [|1|. While it avoids spurious excitations of 
the center of mass, it does not include all spin-orbit pairs of orbitals and thus does not 
obey the Ikeda sum rule for GT strengths. However, this model space (with the choice of an 
appropriate effective interaction) should adequately describe those low-lying states expected 



to be the most important for 2i>f3(3 decay |]T6 . 



We performed the ^^Ge calculation ai (3 = 2.5 MeV~^ with Nt = 60. The effective 
interaction used reproduces the experimental mass splitting of ^^Ge and ^^Se very well: 
20.41 ± 1.3 MeV compared to the experimental splitting of 20.72 MeV (the Coulomb energy 
was calculated following Ref JT^])- Our value for the P~ strength of ^^Ge is B{GTJ) = 
(G^ ■ G) = 12.6±0.3 and we find an energy centroid of 6.3±0.2 MeV, while the experimental 
values are 19.9 and 9.1 MeV respectively JlSl; the discrepancies are almost certainly due to 



not having the gj/2 and /7/2 orbitals in our model space. The missing high-energy strength 
should not affect the low-lying states of ''^As that are important for 2z//5/5 decay. We find a 
vanishingly small f3~^ strength for ^^Se, B{GT^) = —0.05 ±0.03. This strength is identically 
zero in the independent particle model and it is generated only by the smearing of the fermi 
surface due to the interaction. As our (yf-extrapolation takes place over a range of pairing 
strengths, B{GT^) varies by almost an order of magnitude between g = —1.0 and g = 1.0, 
leading to large uncertainties in the physical value. A detailed study of this mass region 
with our effective interaction is planned. 

We performed two independent sets of calculations for both the closure and the exact 
matrix element using the x = 4 and x = oo families of H{g,x)- The extrapolations to 
= 1.0 were linear for both the closure and the exact matrix elements in both cases. Our 
results for the closure matrix elements are —0.36 ±0.34 and 0.08 ±0.17 for x = 4 and x = oo 
respectively. These are to be compared with the truncated shell model calculation of Haxton 
et al. |jl|] (using a different effective interaction) that resulted in a value of 2.56. 

We find consistent exact matrix elements for the x = 4 and x = oo cases (Fig. 3). 
Our results are 0.12 ± 0.07 MeV^^ and 0.14 ± 0.08 MeV^^ respectively (a combined value 
of 0.13 ± 0.05 MeV~^), while the experimental value of this matrix element (using gj^ = 
1.26) is 0.14 ± 0.01 MeV~^ [Q. However, shell model calculations of ordinary /5-decay 
consistently suggest that gA is renormalized to 1.0 in the nuclear medium in which case 
the experimental matrix element is 0.22 ± 0.01 MeV~^. 

There has been no previous shell model calculation of M"^^. Haxton et al. [|I|] obtained 



an estimate in the closure approximation by taking the average energy denominator to be 
the position of the /3~ GT resonance in ^^Ge (9.4 MeV). We find significantly smaller values 
of ^ = Mc/M2^ (-3.0 ± 3.3 MeV and 0.57 ± 1.26 MeV for x = 4 and oo, respectively), 
evidently due to large cancellations among the various terms in (1). 

In this Letter, we have demonstrated an SMMC method to calculate 2i//3/3 decay matrix 
elements. Results for the ^®Ca decay compare well with those from direct diagonalization. 
We have also calculated the matrix element for ^^Ge in a model space significantly larger than 
previous calculations, and obtain a result that is in reasonable agreement with experiment. 
A more detailed description of these calculations will be given elsewhere, and work is in 
progress to calculate the matrix elements of several other, heavier nuclei. 

This work was carried out under grants from the NSF and the DOE. Computational 
cycles were provided by the Concurrent Supercomputing Consortium on the Intel Touchstone 
Delta and the Intel Paragon and on the IBM SP2 at the Maui High Performance Computing 
Center. 
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FIGURES 

FIG. 1. Upper: ln[^(r,0)] for ^^Ca calculated at /? = 2.0 MeV"^ with Nt = 48. The lines are 
best fits. Lower: SMMC and direct diagonalization closure matrix elements for ^^Ca. The SMMC 
points are linearly extrapolated to g = 1.0. 

FIG. 2. Upper: M^''{T,t) for ^^Ca calculated at /? = 2.0 MeV-^ with Nt = 48. The points at 
large T show the asymptotic value (T oo) of the matrix elements; lines are drawn to guide the 
eye. Lower: SMMC exact matrix elements and the direct diagonalization results for *^^Ca. The 
SMMC matrix elements are linearly extrapolated to g = 1.0. 

FIG. 3. SMMC exact matrix elements for ^^Ge calculated using the Hamiltonians H{g, x) with 
X = 4 and x = oo. The lines are linear fits to the points in both cases. The extrapolated values 
and the experimental result of Ref. [5] are shown staggered around g = 1.0 for clarity. 
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